Proyecto - Contingencias de Vida II

Librerías e importaciones

source("cod/setup.R")
source("cod/load_database.R")
source("cod/prob_estim.R")

Factor de degradación para estados CAR

La metodología original empezaba a dar problemas con probabilidades negativas a partir de una edad aproximada de 95 años, por lo que se decidió implementar un factor de reducción desde los 90 años para primero, complementar la probabilidad creciente de muerte y además poder arreglar el problema de probabilidades negativas.

Construir las tablas

source("cod/tablas.R")
Males <- tablas(1)
Females <- tablas(2)

Mejora de mortalidades en el tiempo y mejora de transiciones de empeoramiento

source("cod/degradar_mort.R")
source("cod/calculo_anualidades.R")
edad20 <- degradar_mort(20, 1)

Comprobación de mejoras

calculo_acumulado(20, edad20)
##          Able     Mild Moderate   Severe Profound     Dead
## [1,] 47.37206 6.189601 3.146575 2.593122 3.690927 37.00772
edad20sin_m <- lapply(Males, function(x) as.data.frame(x[21:120,]))
calculo_acumulado(20, edad20sin_m)
##          Able     Mild Moderate   Severe Profound     Dead
## [1,] 42.42812 5.499605 2.619834 1.856993 2.063667 45.53178

Hay una clara diferencia entre mejorías de mortalidades

Cálculo de valores presentes

Se puede realizar varios seguros con los resultados de calculo_vp. Nótese que estamos en edad 20

prueba <- calculo_vp(20, edad20, 0.07, 0.03)

# Seguro de vida normal, 100 millones
(prueba[6]*100e6 )/(12*prueba[1])
## [1] 40566.67
# Seguro de vida con anualidades en caso de Severe o Profound, pagando Mild y Moderate
(prueba[6]*100e6 + 12*(1.5e6*prueba[4] + 3e6*prueba[5]))/(12*(prueba[1]+prueba[2]+prueba[3]))
## [1] 131244.9
# Seguro de vida con anualidades pagando 0.25e6 en aumento de estado
(prueba[6]*100e6 + 12*(0.25e6*prueba[2] +
                         0.5e6*prueba[3] +
                         0.75e6*prueba[4] +
                         1e6*prueba[5]))/(12*prueba[1])
## [1] 111971.3
# Seguro de vida con anualidades pagando 0.5e6 en aumento de estado
(prueba[6]*100e6 + 12*(0.5e6*prueba[2] +
                         1e6*prueba[3] +
                         1.5e6*prueba[4] +
                         2e6*prueba[5]))/(12*prueba[1])
## [1] 183375.8

Cálculo de las primas

source("cod/prima.R")
# primas_h <- sapply(20:70, function(x) prima(calculo_vp(x, degradar_mort(x, 1), 0.05, 0.03)))
# primas_m <- sapply(20:70, function(x) prima(calculo_vp(x, degradar_mort(x, 2), 0.05, 0.03)))
# df_primas <- data.frame(x = 20:70, hombres = primas_h, mujeres = primas_m)
# write_xlsx(df_primas, "res/primas.xlsx")
df_primas <- read_xlsx("res/primas.xlsx")

Portafolio

Generación del portafolio

Se utiliza una normal para centrar las observaciones en una edad de interés

set.seed(70707)
portfolio <- data.frame(edad = round(rnorm(5000, mean = 45, sd = 6.5)),
                         sexo = round(runif(5000, 1, 2))) %>% 
  arrange(., edad, sexo) %>%
  mutate(id = dense_rank(paste(edad, sexo)))
descripcion <- portfolio %>% count(edad, sexo)

Y se genera la lista de probabilidades

lista <- list()
for(i in 1:length(descripcion$edad)) {
  prob_matrices <- degradar_mort(descripcion$edad[i], descripcion$sexo[i])
  
  # Se acumulan las probabilidades previo a la simulación estocástica
  lista[[i]] <- lapply(prob_matrices, function(df) {
    t(apply(df[, 2:7], 1, cumsum))
  })
}

Representación del portafolio

Prima nivelada

prima_n <- function(interes, inflacion){
  primas_p <- sapply(1:length(descripcion$edad),
                   function(x) calculo_vp(descripcion$edad[x],
                                          degradar_mort(descripcion$edad[x],
                                                        descripcion$sexo[x]),
                                          interes, inflacion))
  nivelada <- primas_p %*% descripcion$n
  return(prima(nivelada))
}
t <- proc.time()
(nivelada <- prima_n(0.05, 0.03))
## [1] 266612.7
proc.time()-t
##    user  system elapsed 
##   10.03    0.16   13.05

Análisis de sensibilidad

# names <- paste(as.character(3:7), "%", sep = "")
# tabla <- sapply(3:7/100,  function(x) sapply(1:5/100, function(y) prima_n(x, y)))
# tabla <- data.frame(tabla, row.names = names)
# colnames(tabla) <- names
# write_xlsx(tabla, "res/sensibilidad.xlsx")
tabla <- read_xlsx("res/sensibilidad.xlsx")
tabla
## # A tibble: 5 × 5
##      `3%`    `4%`    `5%`    `6%`    `7%`
##     <dbl>   <dbl>   <dbl>   <dbl>   <dbl>
## 1 264678. 220266. 184793. 156431. 133701.
## 2 320876. 265653. 221434. 186039. 157679.
## 3 389725. 321474. 266613. 222588. 187272.
## 4 473480. 389725. 322062. 267558. 223727.
## 5 574562. 472596. 389725. 322640. 268489.

Modelo estocástico

Proyección de valores presentes

source("cod/proy_prima.R")
# n <- 100000
# set.seed(70707)
# t <- proc.time()
# proy_prima_data <- proy_prima_par(n, 0.05, 0.03)
# proc.time()-t
# raw <- proy_prima_data
# proy_prima_data <- list()
# for(i in 1:n){
#   proy_prima_data[[i]] <- raw[,,i]
# }

Proyección de primas

# primas_est <- sapply(proy_prima_data, function(x) sapply(
#   1:89, function(y) prima(x[y, ])))
# est_mean <- sapply(1:89, function(x) mean(primas_est[x, ]))
# est_perc <- sapply(1:89, function(x) quantile(primas_est[x, ], 0.995))
# datos <- data.frame(x = descripcion$edad, 
#                     sexo = descripcion$sexo,
#                     mean = est_mean,
#                     perc = est_perc) %>%
#   mutate(sexo = case_when(
#     sexo == 1 ~ "Hombres",
#     sexo == 2 ~ "Mujeres"
#   ))
# write_xlsx(datos, "res/primas_est.xlsx")
datos <- read_xlsx("res/primas_est.xlsx")
fig <- ggplot(datos, aes(x = x)) +
  geom_line(aes(y = mean, color = sexo)) +
  labs(title="Esperanza de las primas estocásticas", x="Edad", y="Frecuencia") +
  scale_color_manual(name="Sexo",
                    values=c("Hombres" = "#0073E6", "Mujeres" = "#FF1493"),
                    labels=c("Hombres", "Mujeres"))+
  theme_minimal()
fig %>% ggplotly()
fig <- ggplot(datos, aes(x = x)) +
  geom_line(aes(y = perc, color = sexo)) +
  labs(title="Percentil 99.5 de las primas estocásticas", x="Edad", y="Frecuencia") +
  scale_color_manual(name="Sexo",
                    values=c("Hombres" = "#0073E6", "Mujeres" = "#FF1493"),
                    labels=c("Hombres", "Mujeres"))+
  theme_minimal()
fig %>% ggplotly()

Proyección de prima nivelada

# primas_estniv <- sapply(proy_prima_data, function(x) prima(descripcion$n %*% x))
# write_xlsx(data.frame(primas_estniv), "res/prima_estniv.xlsx")
primas_estniv <- read_xlsx("res/prima_estniv.xlsx")
mean(primas_estniv[[1]])
## [1] 255467
as.double(quantile(primas_estniv[[1]], 0.995))
## [1] 418650.9

Pérdida máxima probable

# pmp <- data.frame(x = descripcion$edad, sexo = descripcion$sexo, pmp = 0)
# for(fila in 1:89){
#   perdida <- c(0,0,0)
#   for(col in 4:6){
#     data <- sapply(proy_prima_data, function(mat) mat[fila, col])
#     perdida[col-3] <- quantile(data, 0.975)
#   }
#   pmp$pmp[fila] <- (perdida[1]*1e6 + perdida[2]*2e6 + perdida[3]*5e6)*12
# }
# write_xlsx(pmp, "res/pmp.xlsx")
pmp <- read_xlsx("res/pmp.xlsx")
as.double((pmp$pmp %*% descripcion$n)/5000)
## [1] 405534528

De referencia, el valor presente de las entradas es de 140e6 por persona.

Preparación para modelar estocásticamente

Variables globales

interes <- 0.07
inflacion <- 0.03
edades <- portfolio$edad
rango <- 120 - min(edades)
v <- (1 + inflacion) / (1 + interes)
v_power <- v^(0:rango)
mujeres <- sum(portfolio$sexo == 2)
hombres <- sum(portfolio$sexo == 1)
sexos <- portfolio$sexo == 1
variables <- c("lista",
                "portfolio",
                "sexos",
                "hombres",
                "mujeres",
                "rango",
                "v_power",
                "proyeccion") 

Funciones

source("cod/proyecciones.R")

Proyeccion grupal de pagos y vivos

Única

set.seed(1)
t <- proc.time()
prueba <- proyeccion()
proc.time()-t
##    user  system elapsed 
##    1.70    0.03    2.25

Paralelizado

# t <- proc.time()
# proyeccion_data <- proyeccion_par(100, cores = 2)
# proc.time()-t

Resumen estocástico

Esperanza

source("cod/resumen_estoc.R")
# t <- proc.time()
# media <- esperanza(proyeccion_data)
# proc.time()-t

Percentil

# t <- proc.time()
# percent.995 <- perc_0_995(proyeccion_data)
# proc.time()-t

Guardar las proyecciones

# write_xlsx(media, "res/media.xlsx")
# write_xlsx(percent.995, "res/percentil.xlsx")

Leer las proyecciones

media <- list(
  read_xlsx("res/media.xlsx", sheet = 1),
  read_xlsx("res/media.xlsx", sheet = 2),
  read_xlsx("res/media.xlsx", sheet = 3),
  read_xlsx("res/media.xlsx", sheet = 4)
)
percent.995 <- list(
  read_xlsx("res/percentil.xlsx", sheet = 1),
  read_xlsx("res/percentil.xlsx", sheet = 2),
  read_xlsx("res/percentil.xlsx", sheet = 3),
  read_xlsx("res/percentil.xlsx", sheet = 4)
)

Gráficos

Ingresos y egresos

sum(media[[3]][[1]])
## [1] 10502128064
sum(media[[3]][[2]])
## [1] 9262581939
sum(media[[4]][[1]])
## [1] 11187690601
sum(media[[4]][[2]])
## [1] 11600273306
sum(percent.995[[3]][[1]])
## [1] 10715653825
sum(percent.995[[3]][[2]])
## [1] 8132304125
sum(percent.995[[4]][[1]])
## [1] 11394721399
sum(percent.995[[4]][[2]])
## [1] 10446976695